library(lme4)
library(MuMIn)
library(multcomp)
data<-read.table(file="SM9 Data.txt", header=TRUE, sep="\t")
data_c <- transform(data,cAge=Age-mean(Age),cWords=Words-mean(Words),cPropositions=Propositions-mean(Propositions),cGood=Good-mean(Good),cBad=Bad-mean(Bad),cSurvival=Survival-mean(Survival),cSocial=Social-mean(Social),cMale_Stereo=Male_Stereo-mean(Male_Stereo),cFemale_Stereo=Female_Stereo-mean(Female_Stereo),cEmotional=Emotional-mean(Emotional))
#testing for influence of moral content on : i) sentences correctly recognised, ii) rated severity of consequences, and iii) percentage of original vignette correctly recalled.
moral_recall<-glmer(Recall~ Morality + (1|Set/Participant/Vignette),family = binomial, nAGQ = 1, data = data_c, control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE, optCtrl = list(maxfun = 100000)))
recognition<-glmer(Correct~ Morality + (1|Set/Participant/Vignette),family = binomial, nAGQ = 1, data = data_c, control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE, optCtrl = list(maxfun = 100000)))
consequences<-glmer(Outcome_Score~ Morality + (1|Set/Participant/Vignette),family = binomial, nAGQ = 1, data = data_c, control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
null_recall<-glmer(Recall ~ (1|Set/Participant/Vignette), data = data_c, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
null_recog<-glmer(Correct ~ (1|Set/Participant/Vignette), data = data_c, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
null_conseq<-glmer(Outcome_Score ~ (1|Set/Participant/Vignette), data = data_c, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
anova(moral_recall,null_recall)
anova(recognition,null_recog)
anova(consequences,null_conseq)
aggregate(Outcome_Score*10~Morality, data=data_c, FUN=mean)
#loading data with exclusions for model fitting
data2<-read.table(file="SM9b Data.txt", header=TRUE, sep="\t")
data_c2 <- transform(data2,cAge=Age-mean(Age),cWords=Words-mean(Words),cPropositions=Propositions-mean(Propositions),cGood=Good-mean(Good),cBad=Bad-mean(Bad),cSurvival=Survival-mean(Survival),cSocial=Social-mean(Social),cMale_Stereo=Male_Stereo-mean(Male_Stereo),cFemale_Stereo=Female_Stereo-mean(Female_Stereo),cEmotional=Emotional-mean(Emotional))
#fitting models
m01<-glmer(Recall~ cAge + Gender + cWords + Propositions + cGood + cBad + cSurvival + cSocial + cMale_Stereo + cFemale_Stereo + Emotion + cEmotional +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m02<-glmer(Recall~ cAge + Gender + Propositions + cGood + cBad + cSurvival + cSocial + cMale_Stereo + cFemale_Stereo + Emotion + cEmotional +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m03<-glmer(Recall~ cAge + Gender + cGood + cBad + cSurvival + cSocial + cMale_Stereo + cFemale_Stereo + Emotion + cEmotional +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m04<-glmer(Recall~ cAge + Gender + cGood + cBad + cSurvival + cSocial + cFemale_Stereo + Emotion + cEmotional +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m05<-glmer(Recall~ cAge + Gender + cGood + cBad + cSurvival + cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m06<-glmer(Recall~ cAge + Gender + cGood + cSurvival + cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m07<-glmer(Recall~ cAge + cGood + cSurvival + cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa",calc.derivs = FALSE, optCtrl = list(maxfun = 100000)))
m08<-glmer(Recall~ cAge + cSurvival + cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m09<-glmer(Recall~ cAge + cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m10<-glmer(Recall~ cAge + cSocial + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m11<-glmer(Recall~ cAge + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m12<-glmer(Recall~ cAge + cSocial + cFemale_Stereo +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
m13<-glmer(Recall~ cSocial + cFemale_Stereo + Emotion +(1|Set/Participant/Vignette),data = data_c2, family = binomial, nAGQ = 1,control=glmerControl(optimizer="bobyqa", calc.derivs = FALSE,optCtrl = list(maxfun = 100000)))
#model selection table
model.sel(m01,m02,m03,m04,m05,m06,m07,m08,m09,m10,m11,m12,m13)
#model averaging
average_model<-model.avg(m08,m09,m10,m11,m13)
#emotion comparison using best model (m09)
m09_mcp<-glht(m09,linfct=mcp(Emotion="Tukey"))
aggregate(Recall*100~Emotion, data=data_c2, FUN=mean)
aggregate(Recall*100~Emotion, data=data_c2, FUN=sd)